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Abstract 


In  ground  surveillance  from  an  airborne  or  space-based  radar  it  is  desirable  to  be 
able  to  detect  small  moving  targets,  such  as  tanks  or  wheeled  vehicles,  within  se¬ 
vere  ground  clutter.  For  operational  moving  target  indication  (MTI)  systems  the 
clutter  filter  coefficients  have  to  be  updated  frequently  due  to  rapidly  changing 
interference  environment.  This  report  examines  the  small  sample  size  perfor¬ 
mance  of  different  fast  fully  adaptive  space-time  processors  (STAP)  and  com¬ 
pares  it  to  the  optimum-detector  performance.  These  recently  proposed  tech¬ 
niques,  named  Matrix  Transformation  based  Projection  (MTP)  and  Lean  Matrix 
Inversion  (LMI),  were  originally  developed  to  provide  fast  man-made  jammer 
suppression  in  large  surface  phased  array  radars  with  many  elements.  For  this 
application  they  have  been  proven  to  operate  with  near-optimum  performance, 
yet  with  a  computational  expense  drastically  reduced  from  that  of  the  optimum 
detector  in  most  practical  cases.  The  investigation  herein  focuses  on  the  perfor¬ 
mance  achieved  when  only  a  few  data  samples  are  available  to  adapt  (update)  the 
clutter  filter  coefficient. 

In  this  report,  the  techniques  are  described  and  a  number  of  simulations  carried 
out.  The  two  applications,  STAP  and  jammer  suppression,  are  similar;  both  are 
required  to  suppress  an  interference  which  is  characterized  by  a  certain  number 
of  dominant  eigenvalues  of  the  sample  space-time  covariance  matrix.  Despite  the 
similarities  the  performance  between  the  two  differs  due  to  the  different  shapes 
of  their  eigenvalue  distribution.  The  LMI  is  shown  to  give  the  best  Signal-to- 
Noise-plus-CIutter  Ratio  (SNCR)  for  a  given  computational  load. 
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Resume 


Ce  rapport  examine  la  performance  de  differents  processeurs  spatio-temporels 
totalement  adaptables  (STAP)  ^  haute  vitesse  pour  les  6chantillons  de  petite 
taille,  et  il  la  compare  k  la  performance  du  detecteur  optimal.  Les  techniques 
proposdes  rdcemment,  appelees  projection  h  base  de  transformation  matricielle 
(MTP)  et  inversion  matricielle  oblique  (LMl),  ont  6l6  mises  au  point  ^  I’origine 
dans  le  but  d’assurer  la  suppression  rapide  du  brouillage  de  source  artificielle 
dans  les  radars  h.  r^seau  en  phase  pr^sentant  une  grande  surface  et  dot^s  de  nom- 
breux  dements.  Pour  cette  application,  il  s’est  av6r6  qu’ils  foumissent  une  per¬ 
formance  quasi-optimale  et  qu’ils  exigent,  dans  la  plupart  des  cas  pratiques,  des 
ressources  de  calcul  consid^rablement  r^duites  par  rapport  au  detecteur  optimal. 
L’dtude  ddcrite  porte  essentiellement  sur  la  performance  obtenue  lorsque  seule- 
ment  quelques  6chantillons  de  donn6es  sont  disponibles  pour  adapter  (mettre  h 
jour)  le  coefficient  de  filtrage  du  clutter. 

Ce  rapport  d^crit  les  techniques  utilisees  et  un  certain  nombre  de  simulations  ef- 
fectu6es.  Les  deux  applications,  STAP  et  la  suppression  de  brouillage,  sont  sem- 
blables,  et  les  deux  sont  n^cessaires  pour  eliminer  un  brouillage  caract6ris€  par 
un  certain  nombre  de  valeurs  propres  dominantes  dans  la  matrice  de  covariance 
spatio-temporelle  de  I’^chantillon.  En  d6pit  des  similitudes,  les  performances 
des  deux  se  distinguent  par  les  formes  differentes  de  leur  repartition  des  valeurs 
propres.  Il  est  d^montre  que  la  LMI  procure  le  meilleur  rapport  signal  sur  bruit 
plus  clutter  (SNCR)  pour  une  charge  de  calcul  donnde. 
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Executive  summary 


Background 

In  ground  surveillance  from  an  airborne  or  space-based  radar  it  is  desirable  to  be 
able  to  indicate  small  moving  targets,  such  as  tanks  or  wheeled  vehicles,  within 
severe  ground  clutter.  The  clutter  echoes  received  by  the  radar  show  a  motion- 
induced  spread  of  the  Doppler  spectrum  which  particularly  masks  out  the  slowly 
moving  targets.  The  problem  is  to  suppress  this  clutter  effectively. 

The  optimum  technique  to  address  this  problem  is  space-time  adaptive  process¬ 
ing  (STAP).  In  STAP  an  array  of  antenna  elements,  N,  as  well  as  a  number  of 
pulses,  M ,  is  used  to  estimate  the  clutter.  The  optimum  processor  is  given  by  the 
inverse  of  the  space-time  covariance  matrix  of  the  clutter.  Since  this  matrix  is  of 
order  NM,  the  full  inversion  is  too  computationally  intensive  to  adapt  the  pro¬ 
cessor  frequently.  By  exploiting  the  fact  that  the  rank  of  the  clutter  covariance 
matrix  is,  in  most  cases,  much  smaller  than  the  order  (due  to  the  redundancy  in 
the  space-time  domain),  it  is  possible  to  develop  techniques  requiring  less  com¬ 
putational  effort. 

Principal  results 

In  this  report,  two  different,  recently  developed,  fast  rank  reduction  techniques, 
the  so  called  subspace  or  projection  methods,  are  applied  to  the  problem  of  clut¬ 
ter  suppression.  These  methods  were  originally  proposed  to  provide  fast  man¬ 
made  jammer  suppression  in  large  surface  phased  array  radars  with  many  ele¬ 
ments.  The  basic  theory  behind  them  is  briefly  reviewed  and  their  performance 
is  studied  via  simulations. 

Significance  of  results 

While  most  investigations  in  the  literature  deal  only  with  the  asymptotic  per¬ 
formance,  i.e.  an  infinite  number  of  data  samples,  studies  of  the  small  sample 
size  behaviour  are  rare.  Finite  sample  size  is  the  practical  case  for  operational 
radar  systems  in  which  the  filtering  coefficients  have  to  be  adapted  frequently 
along  the  range  dimension.  The  techniques  described  in  this  report  are  found  to 
produce  superior  Signal-to-Noise-plus-Clutter  Ratios  (SNCR)  for  a  given  com¬ 
putational  complexity.  The  most  computationally  efficient  algorithm  reduces  the 
requirement  for  high  speed  computation  by  10  —  1000  times  in  many  practical 
cases. 


Christoph  H.  Gierull.  2001.  Analysis  of  the  Small  Sample  Size  Performance  of  Fast 
Fully  Adaptive  STAP  Techniques  for  MTI  Radar.  DREO  TR  2001-079.  Defence  Re¬ 
search  Establishment  Ottawa. 
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Sommaire 


Contexte 

Pour  la  surveillance  air-sol,  il  est  souhaitable  de  pouvoir  d^tecter  de  petites  cibles 
mobiles,  comme  des  blind6s  ou  des  v6hicules  k  roues,  dans  du  clutter  de  sol  in¬ 
tense.  Les  dchos  de  clutter  reus  par  I’a^ronef  montrent  un  6talement  du  spec¬ 
tre  Doppler  caus6  par  le  mouvement,  qui  masque  particuli^rement  les  cibles 
mobiles  lentes.  Le  problbme  consiste  k  supprimer  efficacement  ce  clutter.  La 
technique  optimale  pour  aborder  ce  probleme  est  le  traitement  spatio-temporel 
adaptable  (STAP).  Selon  cette  technique,  un  r^seau  d’616ments  d’antenne,  N, 
ainsi  qu’un  certain  nombre  d’impulsions,  M,  sont  utilises  pour  ^valuer  le  clut¬ 
ter.  Le  processeur  optimal  correspond  k  I’inverse  de  la  matrice  de  covariance 
spatio-temporelle.  Comme  cette  matrice  est  de  I’ordre  N  M,  I’inversion  complete 
n^cessite  trop  de  calculs  pour  I’adaptation  frequente  du  processeur.  En  exploitant 
le  fait  que  le  rang  de  la  matrice  de  covariance  du  clutter  est,  dans  la  plupart  des 
cas,  beaucoup  plus  faible  que  I’ordre  (en  raison  de  la  redondance  dans  le  domaine 
spatio-temporel),  il  est  possible  de  mettre  au  point  des  techniques  exigeant  moins 
d’ efforts  de  calcul. 

Principaux  resultats 

Ce  rapport  d6crit  comment  deux  techniques  diff^rentes  de  reduction  rapide  du 
rang  qui  ont  r6cemment  6t6  mises  au  point,  soit  les  m^thodes  subspatiales  ou  de 
projection,  sont  appliqu^es  au  probleme  de  suppression  du  clutter.  Ces  m6thodes 
ont  6t6  mises  au  point  k  I’origine  dans  le  but  d’assurer  la  suppression  rapide  du 
brouillage  de  source  artificielle  dans  les  radars  k  r6seau  en  phase  pr6sentant  une 
grande  surface  et  dot6s  de  nombreux  61^ments.  Les  principes  sous-jacents  sont 
examines  bribvement  et  les  performances  sont  etudiees  h  I’aide  de  simulations. 

Signification  des  resuitats 

Bien  que  la  plupart  des  Etudes  qui  se  retrouvent  dans  le  documentation  portent 
uniquement  sur  la  performance  asymptomatique,  c.-h-d.  sur  un  nombre  infini 
d’6chantillons  de  donnees,  les  dtudes  comportementales  relatives  aux  echantillons 
de  petite  taille  sont  rares.  Dans  la  pratique,  des  tailles  d’echantillon  finies  sont 
utilis^es  pour  les  syst^mes  radar  op6rationnels  dans  lesquels  les  coefficients  de 
filtrage  doivent  s’ adapter  fr^quemment  dans  la  dimension  de  la  distance.  11 
s’av^re  que  les  techniques  d6crites  dans  ce  rapport  produisent  des  rapports  signal 
sur  bruit  plus  clutter  (SNCR)  sup6rieurs  pour  une  complexity  de  calcul  d6terminye. 
L’algorithme  le  plus  efficient  du  point  de  vue  des  calculs  reduit  de  10  h  1  000  fois 
la  necessity  de  calculs  ^  haute  vitesse,  dans  beaucoup  de  cas  pratiques. 

Christoph  H.  Gierull.  2001 .  Analyse  de  la  performance  de  techniques  spatio-temporelles 

totalement  adaptables  (STAP)  k  haute  vitesse  pour  les  ychantillons  de  petite  taille. 

DREO  TR  2001-079.  Centre  pour  la  Recherche  de  la  Dyfence  Ottawa. 
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1.  Introduction 


In  military  space-based  or  air-to-ground  surveillance  it  is  desirable  to  be  able 
to  detect  small  moving  ground  targets  such  as  tanks  or  wheeled  vehicles  within 
severe  ground  clutter.  The  clutter  echoes  received  by  the  aircraft  show  a  motion- 
induced  spread  of  the  Doppler  spectrum  which  particularly  masks  the  slowly 
moving  targets.  The  problem  is  to  suppress  this  clutter  effectively,  i.e.  down  to 
the  sensor  noise  power  level.  One  very  promising  solution  is  the  use  of  mul¬ 
tiple  antenna  apertures  or  a  single  antenna  aperture  with  multiple  feed  horns. 
The  extra  degrees  of  freedom  offered  by  these  structures  allows  one  to  exploit 
the  information  contained  in  both  the  temporal  (Doppler)  and  angular  (spatial) 
domains  towards  the  goal  of  distinguishing  between  the  clutter  and  the  desired 
target  signals.  The  different  phase  centers  of  the  antenna  provide  the  auxiliary 
spatial  dimension  which  allows  the  use  of  digital  signal  processing  to  suppress 
the  unwanted  clutter  energy.  Good  general  introductions  of  space-time  filtering 
of  data  from  an  array  of  sensors  for  Moving  Target  Indication  (MTI)  are  given  in 
[1],  [2]. 

In  this  report,  two  different  recently  developed  fast  rank  reduction  techniques, 
so  called  subspace  or  projection  methods,  are  applied  to  the  problem  of  clutter 
suppression.  They  were  originally  proposed  to  provide  fast  man-made  jammer 
suppression  in  large  surface  phased  array  radars  with  many  elements.  In  this 
report,  the  basic  theory  behind  them  is  reviewed  and  their  performance  is  studied 
via  simulations.  While  most  investigations  in  the  literature  deal  only  with  the 
asymptotic  performance,  i.e.  an  infinite  number  of  data  samples,  studies  of  the 
small  sample  size  behaviour  are  rare.  Finite  sample  size  is  the  practical  case 
for  an  operational  system  for  which  the  filtering  coefficients  have  to  be  adapted 
frequently  along  the  range  dimension. 

In  chapter  2  the  problem  is  mathematically  formulated  and  the  concept  of  opti¬ 
mum  processing  is  briefly  reviewed.  The  simulation  scenario  used  for  the  analy¬ 
sis  and  comparison  of  the  different  methods  is  described.  In  chapter  3  the  finite 
sample  size  processing  scheme  is  introduced.  This  processing  scheme  is  required 
to  estimate  the  covariance  matrix  which  is  not  known  in  practice.  In  section  3.2, 
the  two  SNCR-optimizing  techniques  SMI/LSMI  are  reviewed,  followed  by  the 
description  of  the  idea  behind  subspace  or  projection  methods.  Since  the  com¬ 
putational  complexity  of  the  best  known  subspace  technique  (the  eigenvector 
decomposition)  is  enormous,  a  much  faster  method,  called  MTP,  is  proposed  in 
the  following  section.  The  shape  of  the  eigenvalue  distribution  of  the  covariance 
matrix  is  exploited  to  develop  a  robust  and  computationally  efficient  algorithm. 
It  is  based  on  a  weighting  of  the  projector  and  is  called  lean  matrix  inversion 
(LMI).  The  computational  complexity  of  the  various  techniques  is  compared  in 
chapter  4. 
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2.  Optimum  Processing 


2.1  Maximization  of  the  SNCR 


Let  a  spatial-time-snapshot  x,  hereafter  called  a  snapshot,  be  modeled  as  the  su¬ 
perposition  ^  where  ^  is  the  sensor  noise  vector  with  mutually 

independent  components  of  power  and  the  time  k\—  ^  is  a  complex 

normal  distributed  random  vector  with  independent  entries  representing  the  am¬ 
plitudes  of  the  clutter  at  time  k.  The  snapshots  are  sampled  at  such  intervals  that 
^  and  are  completely  independent  (see  also  [1]).  In  the  case  of  MTI  the 
snapshots  are  often  gathered  from  the  range  bins  surrounding  a  guard  bin  along 
the  range  direction.  With  these  definitions  in  mind,  the  received  data  vectors 
can  be  understood  as  realizations  of  mutually  independent,  identically  complex 
normal  distributed  random  vectors  X-  The  random  vectors  X  have  mean  0  and 
the  covariance  matrix  R  =  EXX^  =  D  AD^  +  a^I,  where  E  is  the  expectation 
operator  and  I  denotes  the  identity  matrix.  Matrices  are  written  bold,  vectors  are 
underlined  and  (•)^  means  complex  conjugate  transpose. 


As  can  be  demonstrated  using  basic  detection  theory  [3],  the  probability  of  de¬ 
tection  of  a  signal  in  noise  depends  mainly  on  the  SNR.  Therefore,  one  criterion 
for  adjusting  the  filter  vector  w  could  be  the  maximization  of  the  signal-to-noise- 
plus-clutter  ratio  (SNCR)  of  the  filter  output  when  the  processor  is  matched  to 
frequency  ■&.  The  SNCR  of  this  filter  output,  hereafter  denoted  as  p,  is  given  by 


P(vh)  :=  p(w(d))  = 


2E|vv"rf(i3)|2 
'  '  E|w^X|2 


w»Rw 


0) 


where  a  is  the  complex  amplitude  and  d{‘&)  a  given  signal  at  this  frequency.  It  is 
easy  to  verify  that  the  solution  of  (1)  is  given  by 

vVop,  =  aR-*^[(i3)  (2) 

with 

Popt=^W"R‘'^W,  (3) 

where  a  is  an  arbitrary  complex  factor  [4, 5].  The  optimum  processor  can  be  sep¬ 
arated  into  two  steps:  a  decorrelation  of  the  clutter  via  multiplication  with  R“*/^ 
(pre-whitening)  and  a  matched  filtering  with  the  adapted  signal  5  =  R' 


2.2  Simulation  Scenario 

In  order  to  analyze  and  compare  different  suppression  techniques,  the  SNCR 
can  be  used  as  a  measure  of  the  performance  of  the  clutter  filter.  The  SNCR  is 
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normalized  with  respect  to  the  signal-to-noise  ratio  in  the  noise-only  case: 

PW  :=  P(h:(*))  = 

oly  K  vtrj  without  clutter 

i.e.  the  resultant  curve  is  actually  the  SNCR-loss  with  respect  to  the  best  possible 
case  of  no  interference  at  all. 

Fig.  1  shows  the  normalized  SNCR  versus  the  normalized  target  frequency  for 
the  two  special  cases  of,  no  adaptation'  with  w(i&)  =  d{'&)  and,  optimum  pro¬ 
cessing  with  w('d)  =  The  simulation  scenario  used  in  this  report  was 

chosen  as  follows: 

•  Side-looking  geometry 

•  Number  of  channels  N  =\l,'k/2  spacing. 

•  Number  of  Pulses  M  =  12,  DPCA  spacing,  i.e.  (vq/PRF  =  X/2). 

•  Clutter-to-noise  ratio  CNR  =  20dB. 

The  sensor  noise  level  has  been  normalized  to  =  OdB.  One  can  see  that  almost 
no  loss  is  suffered  by  optimum  processing  even  for  target  frequencies  close  to 
the  clutter  notch.  In  contrast,  the  case  of  no  adaptation  equals  the  inverse  of  the 
filter  response  function  and  has  an  average  loss  of  about  —  17dB  over  the  entire 
frequency  range.  If  the  desired  target  signal  perfectly  overlays  the  clutter  =  0, 
the  loss  of  SNCR  is  given  as  —lO*loglO{NM)  —  CNR  =  — 40dB. 

As  wiU  be  shown  later,  an  important  measure  of  the  quality  of  the  clutter  sup¬ 
pression  is  given  by  the  magnitudes  and  shape  of  the  eigenvalue  distribution  of 
the  space-time  covariance  matrix.  Fig.  2  shows  the  corresponding  distribution 
of  the  simulated  covariance  matrix.  Since  the  DPCA-criterion  is  fulfilled  in  this 
side-looking  case,  the  number  of  eigenvalues  larger  than  the  sensor  noise  power 
level  can  be  determined  asN+M-1  =23  (Brennan’s-rule  [1]).  Only  the  first 
50  eigenvalues  are  plotted.  In  order  to  simplify  the  notation  throughout  the  re¬ 
port,  the  order  of  the  space-time  covariance  matrix  will  henceforth  be  denoted 
as  N  and  the  number  of  dominant  eigenvalues  (larger  than  a^)  as  M.  This  nota¬ 
tion  corresponds  to  the  jammer  suppression  case.  Accordingly,  for  STAP,  these 
numbers  have  to  be  substituted  with  NM  and  N+M,  respectively, 

*  When  no  clutter  is  present  and  therefore  no  suppression  is  needed. 
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Figure  1 :  Normalized  signat-to-clutter  plus  noise  ratio  versus  normalized  fre¬ 
quency  for  non-  and  optimum  processing. 


Figure  2 :  Eigenvalue  distribution  of  the  simulated  space-time  covariance 
matrix. 


4 


DREO  TR  2001-079 


3.  Finite  Sampie  Size  Processing 


In  practice  the  exact  (asymptotic)  space-time  covariance  matrix  is  not  known 
and  must  be  estimated  from  a  given  number  of  samples.  With  the  number  of 
snapshots  denoted  by  K,  an  often  used  estimation  is  given  by  the  average  of  the 
dyads 


Assuming  that  the  vectors  Xjt  are  identically  complex  normal-distributed  with 
expectation  zero,  R  is  equivalent  to  the  Maximum-Likelihood-estimator  of  R. 
That  is  known  to  to  be  complex  Wishart-distributed  with  A'-degrees  of  freedom 
and  parameter  matrix  R  [6]. 

A  requirement  for  a  large  number  of  samples  often  proves  untenable  for  practical 
reasons.  For  instance,  the  acquisition  time  of  the  samples  may  be  too  long  or 
the  extra  computational  effort  may  be  wasteful  if  an  acceptable  SNCR  could 
have  been  reached  with  fewer  samples.  For  MTI-applications,  a  large  number  of 
samples  are  often  not  available  since  the  snapshots  are  gathered  from  the  range 
bins  surrounding  the  guard  bin  along  the  range  direction.  Since  the  clutter  is 
non-stationary  due  to  its  range-dependence,  only  a  small  number  of  surrounding 
range  cells  might  be  useful  for  estimating  the  clutter  covariance  matrix  of  the 
cell  in  question. 


3.1  Sample  Matrix  Inversion  (SMI) 


In  an  approach  analogous  to  the  maximization  of  the  SNCR  with  known  covari¬ 
ance  matrix  it  is  possible  to  determine  the  filter  vector  by  maximizing  the  SNCR 
using  only  a  limited  number  of  samples 


\y^d^ 


=  max. 


(6) 


The  optimization  criterion  in  (6)  is  equivalent  to  the  minimization  of  the  clutter 
power  with  limited  sample  number  under  the  constraint  of  a  normalized  filter 
vector  in  target  direction, 

^Rw  =  min  u.  c.  =  1 .  (7) 


The  solution  of  eq.  (7)  can  be  obtained  as  a  multiple  of  the  product  of  the  inverse 
sample  covariance  matrix  and  the  signal  vector 

w  =  aR~‘cl.  (8) 
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This  approach  is  well  known  under  the  name  Sample  Matrix  Inversion  SMI  [7, 
5,  8].  To  assure  the  nonsingularity  of  the  space-time  covariance  matrix,  it  is 
necessary  that  the  number  of  snapshots  at  least  matches  the  order  of  the  matrix 
K>N.  Substituting  the  adapted  filter  (8)  into  eq.  (4)  yields  for  the  SNCR 


P(^)  = 


(9) 


where  the  sensor  noise  level  has  been  set  to  OdB.  In  Fig.  3,  10  realizations  of 
the  SNCR  are  plotted  for  SMI.  The  sample  space-time  covariance  matrix  was 
estimated  with  A”  =  150  snapshots  and  shows  a  relatively  large  loss  («  — 18  dB 
)  over  the  whole  frequency  range. 


Sample  Matrix  Inversion,  K  ==  144 


Figured:  Normalized  signal-to-clutter  plus  noise  ratio  versus  normalized  fre¬ 
quency  for  SMI. 


Here,  p(w)  is  also  a  random  variable  depending  on  the  statistics  of  R~*.  For  the 
related  application  of  man  made  jammer  suppression,  theoretical  performance 
analysis  has  been  done  [5].  To  study  the  performance  of  the  SMI,  Reed  et.  al.  [5] 
calculated  the  pdf  of  the  normalized  version 

p(w)  =  p(w)/rf(i>)^R-^rf(i3). 

The  classical  result  they  found  was  that  the  normalized  SNCR  for  SMI  is  Beta- 
distributed  with  the  parameters  {2{N  —  l),2{K  —  N-^2))  and  an  expectation 

Ep  =  (A:-iv-b2)/(/ir+ 1). 
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One  of  the  main  disadvantages  of  the  SMI  is  the  large  number  of  samples  re¬ 
quired  to  adapt  the  weight  vector  and  suppress  the  interference  effectively.  For 
example,  the  average  number  of  snapshots  to  reach  the  3  dB  loss  of  normalized 
SNCR  Ep  =  0.5  is  given  by  K^dB  =  2N-3^2N  (see  also  [7]). 

3.2  Loaded  Sample  Matrix  Inversion  (LSMI) 

3.2.1  Motivation 

In  order  to  improve  the  slow  convergence  of  the  SMI  shown  in  Fig.  3 
one  can  try  to  introduce  an  additional  constraint  to  the  optimization  task 
in  eq.  (7).  For  instance,  y^w  =  const,  can  be  used;  this  also  solves 
another  disadvantage  of  the  SMI,  namely  strong  fluctuating  sidelobes  of 
the  filter  output.  Fig.  4  shows  several  realizations  of  the  adapted  filter 
output  for  the  SMI 

a(i3)  =  \wi-&Qfdi'd)\'^  =  (10) 

The  filter  is  matched  to  a  desired  signal  at  =  0.25  which  means  the 
target  is  located  completely  outside  of  the  clutter  band.  The  bold  curve 
corresponds  to  the  case  where  the  covariance  matrix  is  perfectly  known. 
Obviously,  the  required  clutter  notch  does  exist  but  the  side  lobes  of 
the  filter  output  are  now  strongly  fluctuating  and  greatly  increased  over 
the  entire  frequency  range.  Low  side  lobes  of  the  filter  response  are 
crucial  if  the  probability  of  false  alarm  detections  not  to  increase.  The 
criterion  for  the  minimization  of  the  clutter  power  under  the  constraints 
of  a  normalized  filter  on  one  hand,  and  a  preserved  constant  average 
sidelobe  level  of  the  adapted  filter  response  on  the  other  hand,  yields 

p(w)  =  vt^Rw  =  min  u.  c.  v^^=l,  J  \v^ d{^(C}^dK  =  c.  (11) 
The  integral  in  eq.  (1 1) 

J  \}^d{K)\^dK  =  d{K)d{K)^d\^  w  =  y^Cw 

is  performed  over  the  surface  of  the  unity-sphere.  For  the  case  of  plane- 
antennae  with  equidistant  spacing  of  it  has  been  shown  that  C  con¬ 
verges  quickly  towards  the  identity  matrix.  For  a  linear-antenna  with 
equidistant  spatial  spacing  |  the  identity  C  =  I  holds  exactly  [9].  Us¬ 
ing  this  simplification  and  omitting  the  indices,  ifl,  in  the  following  for 
simplification,  eq.  (1 1)  can  be  expressed  as 

p(w)  =  =  min  u.  c.  v^d^l,  y^yv  =  c  (12) 
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or  as  the  Lagrange-function 

pfejW)  Y)  =  + niy^d  - 1 )  +  y(vh^vv  -  c) , 

respectively. 


Normalized  frequency 


Figure  4 :  Adapted  antenna  pattern  with  SMI,  A"  =  150. 


The  solution  results  in 


w  =  a(R+Yl)  ^d,  (13) 

where  y  must  be  determined  under  the  second  constraint  in  eq.  (12). 
Unfortunately,  no  simple  analytic  relation  between  y  and  c  is  known. 
Therefore,  the  additive  loading  of  the  diagonal  of  the  sample  covariance 
matrix  with  a  suitable  chosen  constant  5,  i.e.  w  =  a  (]^+ 61)  d,  will 
be  taken  as  the  solution  [10,  8,  1 1] .  A  definition  of  a  reasonable  range 
of  values  is  given  in  [12,  13]. 

Fig.  5  shows  the  SNCR  for  the  filter  vector  calculated  with  LSMI- 
technique.  The  number  of  snapshots  was  chosen  to  be  =  M  =  24 
and  the  loading  factor  is  set  to  three  times  the  sensor  noise  power  level, 
i.e.  8  =  3a^. 
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Adapted  antenna  patlern[dB]  ^  Normalized  SNCR  [dB] 


Loaded  Sample  Matrix  Inversion,  K  =  24 


-0.5  -0.4  -0.3  -0.2  -0.1  0  0.1  0.2  0.3  0.4  0.5 


Normalized  Frequency 


re  5 : 


Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for 


LSMI. 


Normalized  frequency 


Figure  6 :  Adapted  antenna  pattern  for  LSMI,  =  150,  5  =  3. 
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A  significant  performance  improvement  over  that  of  the  original  SMI 
can  be  recognized  over  the  entire  frequency  range.  The  approximate 
SNCR  loss  is  only  —2.5  dB,  although  the  number  of  samples  chosen 
was  six  times  less  than  for  Fig.  3.  The  improvement  of  the  SNCR  is 
linked  to  the  desired  reduction  of  fluctuation  in  the  average  sidelobe 
level  of  the  adapted  matched  filter  response,  as  shown  in  Fig.  6,  [14]. 


3.2.2  Consequences 

To  study  the  impact  of  the  diagonal  loading  on  the  adapted  matched 
filter  response,  i.e.  the  gain  of  SNCR,  let  the  sample  covariance  ma¬ 
trix  R  be  decomposed  into  its  real  eigenvalues  (with  the  smallest 
eigenvalue  ^in  >  0)  and  corresponding  orthonormal  eigenvectors  u„ 


N  N 

R  =  with  Y,  =  I- 

n=\  n=\ 


The  inverse  can  be  written  as 


iv  1  \  f  N  \  1 

R-'=  = 

Amin  J 


(14) 


Substituting  eq.  (14)  into  eq.  (10)  yields  the  adapted  filter  response 
without  loading: 


a(i&)  = 


1 


For  a  given  target  signal  orthogonal  to  the  clutter  subspace,  with  rf(i&)  ± 
{«i )  •  •  • )  “a/}’  the  sidelobe  regions  of  the  clutter,  it  holds  that: 


a(i3)  = 


N 


n-dofm-  E 

rt=M+l 


iin 


K 


d{‘&o)^M„y^d{'&) 


When  the  noise  eigenvalues  are  of  different  magnitude,  and  in  partic¬ 
ular  when  the  sample  sizes  are  small,  the  error  term  ^1  —  for 
n  e  {M+  can  be  large.  On  the  basis  of  this  notation,  the  effect 
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of  the  diagonal  loading,  i.e.  =  X.„  +  8,  becomes  clear: 


a{'&)  = 


n=M+l  \  ^ 


(^min  "b  S)  N 
+  5  ) 


’d(po)^Un^d.{‘^) 


1 


(^roin  +  5)^ 
The  error  term  after  diagonal  loading 


g 2  _  ^min  ~b  ^ ^roin 

A,7j  +  8  X,/|  +  8 

becomes  smaller  due  to  the  compression  of  the  noise  eigenvalues,  if  8 
is  larger  than  the  largest  noise  eigenvalue  X.m+i.  On  the  other  hand,  to 
preserve  the  clutter  notch,  8  must  be  chosen  to  be  much  smaller  than 
the  smallest  dominant  clutter  eigenvalue,  i.e.  8  «C  A-m-  For  the  case  of 
adaptive  jammer  suppression,  a  less  known  analytical  expression  for  the 
pdf  of  the  normalized  SNCR  resulting  from  LSMI  can  be  found  in  [12]. 
Herein  it  is  shown  that  if  >  M ,  and  some  reasonable  conditions  are 
imposed  on  the  diagonal  load  y  the  normalized  SNCR  is  also  Beta- 
distributed  but  with  parameters  {2{K  —  M+1),  2M) 

/p(p)  =  y{K-M+  l)r{M)^  0<P<1  (15) 

and  expectation 

Ep  =  (^-M+l)/(^+l).  (16) 


It  can  be  seen  that  the  expected  3  dB  loss  of  SNCR  is  now  reached  for 
only  KidB  =  2A/  —  1  snapshots,  a  number  remarkably  lower  than  that  re¬ 
quired  with  SMI.  Therefore,  the  statistical  properties  are  independent  of 
the  number  of  channels,  and  the  number  of  snapshots  needed  to  suppress 
the  interference  effectively  is  now  proportional  to  M  (see  also  [15, 13]). 

To  study  the  statistical  properties  of  the  improvement  factor  resulting 
from  the  application  of  clutter  suppression  in  more  detail,  the  empir¬ 
ical  probability  density  function  has  been  computed  at  one  particular 
frequency.  In  Fig.  7  the  histogram  of  the  norm^zed  SNCR  generated 
from  100000  trials  is  plotted  along  with  the  theoretical  pdf  in  eq.  15 
for  M  —\2.  The  discrepancy  between  the  two  curves  is  due  to  the  dif¬ 
ferent  shapes  of  the  eigenvalue  distributions.  While  the  eigenvalues  in 

<  Y<  X,Af,  where  Xm  denotes  the  smallest  interference  eigenvalue 
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the  man-made  jammer  case  resemble  step-functions  (sharp  edge),  the 
eigenvalues  in  Fig.  2  decrease  gradually  down  to  the  level  of  the  sen¬ 
sor  noise.  This  fact  is  discussed  in  more  detail  in  the  following  section. 
It  is  noteworthy  to  point  out  that  the  variance  (variation)  of  the  SNCR 
for  clutter  suppression  is  remarkably  smaller  than  in  the  case  of  jammer 
suppression. 


Loaded  Sample  Matrix  Inversion,  K  =  24 


Figure  7 :  Histogram  of  SNCR  for  LSMI  at  normalized  frequency  0.25,  calculated 
over  100000  trials  (dashed).  Theoretical  pdf  (15)  for  M=  12  (solid). 


3.3  Subspace  Techniques 

Using  the  low  rank  approximation  of  the  space-time  covariance  matrix,  subspace 
techniques  can  be  exploited  for  complete  rejection  of  the  clutter  in  terms  of  inde¬ 
pendency  of  clutter  power.  The  basic  idea  of  the  so  called  projection  methods  is 
the  separation  of  the  overall  space  of  observation  into  an  interference  and  a  noise 
part,  followed  by  a  projection  of  the  array  outputs  into  the  clutter-free  subspace 
to  remove  the  interference. 

The  adaptive  spatial  interference  cancellation  problem  is  in  principle  closely  re¬ 
lated  to  that  of  superresolution  angle  estimation.  The  application  of  a  projec¬ 
tion  into  the  space  complementary  to  the  interference  or  jammer  subspace,  often 
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called  the  (sensor)  noise  subspace,  generates  deep  notches  in  the  antenna  pattern 
in  the  directions  of  the  jammers.  The  detection  of  these  notches  leads  directly 
to  estimations  of  the  unknown  locations.  Thus,  all  subspace  techniques  can  be 
used  for  both  problems,  although  their  statistical  properties  can,  of  course,  be 
very  different.  An  overview  of  subspace  methods  concerning  both  applications 
can  be  found  in  [16, 17]. 

Starting  from  the  model  in  chapter  2,  a  clutter  data  vector  can  be  written  as 


x  =  T>b  +  n 


the  superposition  of  the  phase  vectors  D  =  [rfj , . . . ,  weighted  with  the  clutter 
amplitudes  b,  and  the  sensor  noise  n.  Via  a  projection  orthogonal  to  D 

X  =  P^x, 


(i.  e.  the  M-dimensional  subspace  of  the  (C  which  is  spanned  by  the  column 
vectors  rfj , . . .  ,^+a/)  it  is  possible  to  suppress  the  clutter  exactly  down  to  zero. 
Then,  the  filtering 

y  =  =  d^P^D^ 

gives  the  filter  vector  as 


vv  =  P^d. 


For  the  special  case  of  an  indefinite  clutter-to-noise  ratio  CNR  this  orthogo¬ 
nal  projection  is  equivalent  to  the  SNCR  optimization  w  =  Using  the 

matrix  inversion  lemma^,  the  inverse  of  the  space-time  covariance  matrix  R  = 
DBD^  -I-  <3^1  can  be  written  as 

R-’  =  ^(l-D(D^D-bB-*)''D^), 

where  B  =  B/a^.  Assuming  CNR  — > it  follows  B“*  — ^  0,  so  that 

a^R-i ->P^  =  I-D(D^D)“*D".  (17) 


Since  the  phase  matrix  D  is  generally  unknown,  the  projection  matrix  has  to  be 
estimated.  This  can,  for  instance,  be  achieved  by  a  subspace  fitting  of  the  random 

^  NxM  ^  N  ,  r  , 

data  vectors.  A  M-dimensional  subspace  S  e  G  of  the  G  must  be  found  so 
as  to  best  approximate  the  data  vectors  in  a  Mean  Square  Error  MSE  sense: 

E||X-PsX|P  =  min.  (18) 

3  (A  -t-  BCB«)  ~  ‘  =  A- '  -  A- '  B  (B" A-  •  B  -t-  C-  ‘ ) " '  B« A"  ‘ 
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Pg  denotes  a  projection  matrix  into  the  subspace  which  is  spanned  by  the  columns 
of  S  =  [5] , . . .  The  minimization  in  eq.  (18)  is  equivalent  to  the  maximiza¬ 
tion  of 

E||Ps^ll^=  EX^Ps2C=  tr{PsEXX^}  =  tr{PsR}  =  niax. 

If  the  orthonormality  of  the  vectors  5] , is  also  demanded,  the  projection 
can  be  written  as  Pg  =  SS^.  Hence,  the  maximization  turns  out  to  be  an  opti¬ 
mization  under  constraint: 

tr{S^RS}=max  u.  c.  S^S  =  I.  (19) 

Substituting  the  eigenvalue  decomposition  R  =  UAU^  into  eq.  (19)  results  in 

tr{S^UAU^S}  =  tr{Q"AQ}  =max  u.  c.  Q^Q  =  I.  (20) 

The  solution  of  eq.  (20)  is  shown  in  the  appendix  to  be  the  subspace  Q  which  is 
spanned  by  the  eigenvectors  corresponding  to  the  M  dominant  eigenvalues,  or  an 
arbitrary  unitary  transformation  of  it.  Denoting  the  eigenvectors  corresponding 
to  the  M  dominant  eigenvalues  as  S,  the  projection  matrix  can  be  written  as 

P^  =  I-SS^ 

and  the  filter  vector  as 

w  =  oP^d  =  a{d-SS^d).  (21) 


Unlike  the  LSMI,  the  dimension  of  the  clutter  subspace  M  needs  to  be  known 
in  order  to  calculate  Pg  in  eq.  (21);  it  has  to  be  estimated  as  well.  Several 
different  approaches  exist  in  literature.  The  extended  sphericity -test  [18],  tests 
the  equality  of  the  last  N  —  M  smallest  eigenvalues  of  the  sample  covariance 
matrix.  The,  so-called,  information  theory  based  algorithms,  such  as  AIC-  [19] 
and  MDL-criterion  [20]  are  different  approximations  of  the  density  function  of 
the  extended  sphericity-test  statistic.  The  white  noise  test  for  active  systems  [21] 
is  based  mainly  on  the  comparison  of  the  estimated  noise  power  level  with  a  a- 
priori  known  level.  The  analysis  of  a  suitable  dimension  of  the  clutter  subspace 
is  a  major  task  of  this  report. 

3.3.1  Eigenvector  Projection  (EVP) 

Since  the  asymptotic  covariance  matrix  R  is  not  known  in  practice,  it 
is  possible  to  use  the  eigenvectors  of  the  sample  covariance  matrix  to 
determine  the  filter  vector  w 

R  =  UAU^  =  SAsS^  +  NN". 
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Using  eq.  (5),  one  can  show  that  the  desired  subspace  S  G  C  is  the 
one  which  best  approximates  the  sample  vectors  ^  (A;  =  1 , . . . ,  i?")  in  the 
meaning  of  least  square  error  (compare  also  eqs.  (18)  to  (20)): 

K 

]Clb-Ps^|P  =  min 

^=1 


or 


K 


I 


llPsalP 


“■js"  £ife£E's|  =  i:tr{S''RS} 

U.C.  S^S  =  Im, 


max 


(22) 


respectively.  This  technique  is  completely  analogous  to  the  Multiple 
Slgnla  Classification  MUSIC  algorithm  for  high  resolution  spectral  es¬ 
timation  [22].  It  has  been  shown  analytically  that  the  suppression  per¬ 
formance  of  the  so  called  Eigenvector  Projection  (EVP)  is  equivalent 
to  that  of  LSMI  [23].  To  be  consistent,  the  eigenvectors  of  the  sam¬ 
ple  covariance  matrix  are  denoted  by  S  and  the  corresponding  estimated 
weight  vector  by  w('fi)  =  (l  —  SS^)  rf('fi).  Inserting  this  weight  into  eq. 
(1),  we  get  for  the  SNCR 


^(i3)^P^RP^d(d)  N  ' 


(23) 


The  behaviour  of  p  is  shown  in  Fig.  8  and  Fig.  10  for  a  clutter  subspace 
dimension  of  12  and  24,  respectively.  One  can  see  that  the  average  loss 
of  SNCR  in  both  cases  is  about  (—2)  to  (—3)  dB.  To  evaluate  the  statisti¬ 
cal  properties  of  the  SNCR  in  more  detail,  identical  histograms  to  those 
used  with  LSMI  have  been  calculated  and  displayed  in  Fig.  9  and  11. 
Two  remarkable  results  are  observed. 


In  regards  to  asymptotic  losses,  it  can  be  seen  that  even  a  small  dimen¬ 
sion  of  the  clutter  subspace  results  in  near-optimum  suppression.  The 
main  difference  is  some  minor  ripples;  the  reason  for  these  ripples  lies 
in  the  special  shape  of  the  eigenvalue  distribution.  The  sum  of  all  eigen¬ 
values  is  identical  to  the  total  clutter  power.  After  removing  the  power 
contained  in  the  first  12  eigenvalues  the  power  in  the  remaining  eigen¬ 
values  is  only  a  small  fraction  of  the  entire  clutter  power,  even  if  they  are 
still  larger  than  the  sensor  noise  floor.  Further  increasing  the  subspace 
dimension  to  24  (Fig.  10)  eliminates  the  remaining  power  and  causes 
the  ripples  to  disappear. 
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Normalized  SNCR  [dB] 


In  contrast,  when  a  limited  number  of  samples  is  used  to  estimate  the 
covariance  matrix,  surprising  behaviour  is  observed,  namely  the  smaller 
subspace  dimension  produces  slightly  better  results,  i.e.  the  loss  appears 
to  be  slightly  less.  This  impression  is  confirmed  by  the  two  histograms 
in  Fig.  9  and  11,  where  the  demonstrated  expectation  is  slightly  smaller. 
The  reason  for  this  behaviour  is  not  clear  but  it  is  believed  that  both  the 
magnitude  of  the  eigenvalues  and  the  orientation  of  the  corresponding 
clutter  and  noise  subspaces  with  respect  to  each  other,  play  an  important 
role  (compare  also  section  3.3.4). 


Eigenvector  Projection,  SS  =  12,  K  *  24 
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Figure  8 :  Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for  EVP, 
dimension  of  subspace  12. 
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Histogram 
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Figure  11 :  Histogram  of  SNCR  for  EVP  at  normalized  frequency  0.25  calculated 
over  10000  trials,  dimension  of  subspace  24. 


3.3.2  Projections  without  Eigenanaiysis 

The  estimation  of  the  complete  covariance  matrix  and  the  following 
eigendecomposition  are  very  intensive  and  time  consuming  computa¬ 
tional  processes.  Many  different  techniques  have  been  proposed  in  the 
literature,  with  the  goal  of  reducing  the  computational  complexity  of 
subspace  methods  for  use  in  real  time  applications  in  large  antenna  ar¬ 
rays,  see  e.g.  [24,  25,  26].  These  iterative  methods,  which  employ  ex¬ 
ponential  forgetting  factor  to  partially  update  the  eigenvectors,  have  the 
drawback  of  possessing  a  lower  convergence  rate  than  MUSIC,  i.e.  they 
need  a  large  number  of  samples  to  perform  sufficiently  well.  In  con¬ 
trast,  subspace  methods  which  work  with  successive  short  time  frames 
(i.e.  on  a  small  block  of  data)  show  better  tracking  performance  [17]. 

Since  every  arbitrary  complete  basis  of  the  clutter  subspace  can  be  used 
to  form  the  projection  matrix,  it  is  possible  to  develop  faster  subspace 
methods  which  are  not  based  on  eigenanaiysis.  In  principle,  the  simplest 
way  to  estimate  the  subspace,  as  suggested  by  Hung  and  Turner  [27],  is 
to  take  the  data  vectors  directly  as  a  basis  for  the  clutter  subspace.  This 
method,  called  HTP,  gives  only  a  rough  estimate  of  the  subspace  and 
yields  sufficient  results  only  for  high  CNR.  The  theoretical  statistical 


Eigenvector  Projection,  SS  =  24,  K  «  24 
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analysis  of  the  HTP,  when  used  for  man-made  jammer  suppression,  can 
be  found  in  [13]. 

The  noneigenvector  estimation  of  a  basis  of  the  signal  subspace  based 
on  the  sample  space-time  covariance  matrix,  can  be  conceptualized  as 
an  extension  of  the  HTP  to  the  case  of  an  infinite  amount  of  samples. 
One  class  of  subspace  estimation  methods,  called  Yeh-Brandwood-  Pro¬ 
jections  (YBP),  is  based  on  the  idea  of  a  sufficient  data  reduction  pre¬ 
processing,  namely  the  estimation  of  a  basis  of  the  clutter  subspace 
based  only  on  a  part  of  the  sample  covariance  matrix.  Since  the  mea¬ 
sured  covariance  matrix  always  contains  receiver  noise,  the  estimated 
subspace  is  disturbed. 

Different  techniques  have  been  proposed  in  the  literature  to  reduce  the 
negative  influence  of  the  sensor  noise  e.g.  [28, 29].  In  [30]  it  was  shown 
that  the  influence  of  the  additive  noise  can  be  interpreted  as  a  leakage  of 
the  signal  power  into  the  noise  subspace.  Instead  of  recollecting  the  lost 
jammer  power,  it  is  desirable  to  pre-apply  an  adaptive  mechanism  which 
protects  the  method  before  or  during  execution  of  the  subspace  estima¬ 
tion.  A  mathematical  technique  to  achieve  this  property  is  known  as  the 
rank-revealing  QR-factorization  [31]  of  the  sample  covariance  matrix, 
also  called  QR-decomposition  with  permutation  (PQR)  [32].  This  tech¬ 
nique  was  first  used  in  this  context  by  Hsieh  et.  al.  [33]  and  Prasad  et.al 
[34]. 

To  avoid  the  time  consuming  recursive  permutation,  Reilly  [35]  suggests 
a  fixed  pre-selection  of  A/  columns  of  the  covariance  matrix  using  some 
a-priori  knowledge  about  their  structure.  He  called  this  method  reshuf¬ 
fled  QR  decomposition,  or  RQR.  In  [36],  a  simple,  efficient  extension 
of  the  RQR,  called  Matrix  Transformation  based  Projection  (MTP),  has 
been  presented.  It  is  based  on  a  generalised  representation  of  the  class  of 
YBP.  The  performance  of  the  MTP,  when  applied  to  the  closely  related 
problem  of  superresolution  angle  estimation  with  small  sample  size,  has 
been  investigated  in  [37]. 

3.3.3  Projections  Based  on  Covariance  Matrix  Transformation 
(MTP) 

This  section  deals  with  the  performance  of  MTP  when  applied  to  the 
suppression  of  clutter  for  airborne  MTI.  Let  us  first  consider  a  general 
dimension  reducing  transformation  of  the  covariance  matrix 

R  =  RT  =  (DAD"  -K  0^1)  T,  (24) 
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NxM 

where  T  is  an  arbitrary  transformation  matrix  of  the  ®  with  rank{T } 

=  M.  It  can  be  directly  seen  from  eq.  (24)  that  a  substitution  of  T  by  the 
first  M  columns  of  the  identity  matrix  is  identical  with  the  classical  YBP. 
Moreover,  if  we  replace  the  transformation  by  the  permutation  or  pre¬ 
selection  matrices  mentioned  in  the  previous  section,  the  PQR  and  RQR 
methods  are  also  included  as  special  cases  of  this  general  expression. 
Using  R  =  XX^T  to  build  up  the  projection  matrix 


and  inserting  the  resulting  weight  vector  w(0)  =  P^£[('d)  into  eq.  (1), 
we  get  for  the  SNCR 


.  a(«)''P4rf(«)P  <,2 

p(R,«)  -  H  ■ 

R  R 


(25) 


3.3.3. 1  Optimum  Transformation  for  Known  Covariance  Matrices 

Eq.  (24)  raises  the  question  of  whether  a  constant  (non  adap¬ 
tive)  transformation  T  exists  for  which  an  efficient  subspace 
estimation  can  be  achieved,  assuming  the  covariance  struc¬ 
ture,  including  the  sensor  noise  term,  is  given.  That  means, 
one  is  looking  for  a  transformation  of  the  covariance  ma¬ 
trix  for  which  lin{RT}  =  lin{D},  i.e.  RT  must  be  a  linear 

combination  of  the  phase  matrix  R  =  RT  :=  DG  which  is 

. ,  ^  ^MxM 

not  at  all  affected  by  the  sensor  noise,  with  G  €  G  and 
rank{G}  =  M.  It  is  easy  to  verify  that  an  optimum  transfor¬ 
mation  in  this  sense  is  given  by  an  arbitrary  linear  combina- 

MxM 

tion  of  the  phase  matrix  itself:  T  =  DF  with  F  €  C  and 
rank{F}  =  M.  Since 

R  =  RT  =  RDF  =  DAD^DF  -b  a^DF  =  DG, 
where 


rr  o  MXM 

G  =  AD^DF  +  a^FeG 

R  represents  a  basis  of  the  interference  subspace  which  is  not 
interrupted  by  the  noise  even  in  the  case  of  fully  correlated 
interference  sources,  i.e.  when  rank{  A}  <  M. 
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3.3.3.2  Suboptimum  Transformation  for  Estimated  Covariance  Ma¬ 
trices 

As  the  phase  matrix  of  the  clutter  is  supposed  to  be  unknown, 
one  might  follow  the  considerations  of  the  HTP  and  take  M 
snapshots  vectors  Xq  =  [xj, . . .  ,^]  as  a  rough  estimation  of 

the  interference  subspace  for  the  transformation  R  =  RT  = 
XX^Xq.  The  result  of  this  approach  shows  an  unexpected 
effect.  Although  the  transformation  matrix  T  =  Xq,  i.e.  the 
M  snapshots  are  highly  fluctuating  and  therefore  the  spanned 
subspace  lin{Xo}  strongly  deviates  from  the  exact  one'*  lin{D}, 

it  seems  that  the  subspace  lin{R}  (spanned  by  R  =  D  AD^Xq 
+a^Xo)  has  effectively  no  loss  of  SNCR  compared  to  lin{D}. 
To  explain  this  phenomenon,  we  will  establish  a  relationship 
between  the  MTP  and  the  Power-method  [32]  for  recursive 
evaluation  of  the  eigenvectors  corresponding  to  the  M  domi¬ 
nant  eigenvalues: 

Yi  =  RY/_i  =  R^Yo 

Y/=  orth{Y/}  /=1,...,L,  (26) 

where  orth  {$}  denotes  the  orthonormalization  of  the  columns 
included  in 

One  can  see  in  eq.  (26),  that  depending  on  the  initial  matrix 
Yq  used  in  the  first  step  of  iteration  L  =  1,  either  a  YBP  with 
Yq  =  i,  a  PQR  or  a  RQR  will  be  performed.  To  analyse  the 
convergence  of  the  algorithm,  we  will  replace  the  covariance 
matrix  in  (26)  by  its  eigen-factorization.  Further,  let  the  initial 
matrix  Yq  be  an  arbitrary  linear  combination  of  the  eigenvec¬ 
tors 


Yo  =  UA  =  SAs-bNAA^ 

r.T  .xiT  .  MxM  ,  .  (N-M)xM 

with  A  =  [Aj,A^]  ,  As  €  G  and  €  G 
Assuming  linear  independence  of  the  eigenvectors  and  the 
columns  of  Yq,  it  directly  follows  that  the  matrix  of  coeffi¬ 
cients  A  =  U^Yq  has  the  rank{A}  =  M  and,  particularly, 
that  rank{As}  =  rank{S^Yo}  =  M.  Inserting  these  terms 
into  the  iteration  (26)  yields 

Yi  =  R^Yo  =  S  A^As  +  NA^A/v 

'‘Compare,  for  instance.  Fig.  1  in  [13] 
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and  due  to  the  nonsingularity  of  As 


[s  +  NAjvAyvA^'A^' 
[s+0(^')]AiA5 


A^A^ 
/=!,.. .,L 


(27) 


with 


q  =  max{XA/+i  /X,i , . . . ,  A,m+i  •  •  •  >  ^yv/^i  j 

Am- 


The  convergence  of  the  iteration  is  ensured  by  the  succes¬ 
sive  orthonormalization  of  Y/.  If  we  take  into  account  that 
the  smallest  clutter  eigenvalue  Xm  is  proportional  to  the  av¬ 
erage  of  the  interference  power  of  the  clutter,  and  that  X,a/+i 
is  approximately  equal  to  the  sensor  noise  power  then  the 
order  of  the  error  quotient  is  proportional  to  (1/CNR)^. 
This  can  be  negligible  even  for  L  =  1,  particularly  in  case 
of  strong  clutter.  Thus,  it  can  be  concluded  that  the  perfor¬ 
mance  of  the  MTP  and  YBP  depend  not  only  on  the  CNR, 
but  also  on  the  condition  of  the  coefficient  matrix  A5.  This 
statement  can  easily  be  confirmed  if  a  pure  noise  matrix  for 
transformation  is  used.  In  that  case,  the  coefficient  matrix 
As  =  S^T  is  given  as  a  random  linear  combination  of  the 
clutter  eigenvectors,  which  is  generally  well  conditioned.  Al¬ 
though  this  kind  of  transformation  has  no  relation  to  the  inter¬ 
ference  subspace,  the  given  results  shown  no  difference  (see 
also  [15]).  Under  MTP,  the  number  of  complex  operations 
to  estimate  the  interference  subspace  R  =  RT  =  X  (X^T) 
is  of  the  order  0{2KN).  This  is  double  that  of  the  YBP  or 
RQR-method  which  achieve  greater  computational  efficiency 
through  the  special  structure  of  their  transformations,  e.g.  us¬ 
ing  the  identity  or  any  other  permutation  matrix  for  transfor¬ 
mation  requires  no  explicit  matrix  multiplication  rather  only 
a  selection  of  corresponding  columns. 

3.3.3.3  Suboptimum  Transformation  with  Smaii  Computationai 
Compiexity 

In  this  section,  we  use  a  physical  interpretation  of  the  different 
special  cases  of  the  YBP  introduced  above,  to  develop  a  trans¬ 
formation  matrix  with  very  small  computational  effort.  In 
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other  words,  we  are  looking  for  a  transformation  matrix  which 
is  sparse  but  still  guarantees  a  good  condition  of  As  =  S^T. 

The  estimation  of  a  part  of  the  sample  space-time  covariance 
matrix  corresponds  to  the  average  of  the  correlation  of  a  M- 
element  subarray  output  with  the  entire  array.  For  example,  a 
subarray  built  with  the  first  M  elements  (YBP)  shows  a  wide 
beamwidth  due  to  its  strongly  decreased  aperture.  A  pair  of 
closely  located  clutter  points  which  are  still  conventionally 
resolvable  by  the  entire  aperture  Am  >  beamwidth  ~  l/AT, 
now  appear  so  close  together  that  the  corresponding  subar¬ 
ray  phase  matrix  D^T  in  eq.  (24)  will  be  almost  singular. 
It  now  becomes  clear  that  a  fixed  pre-selection  scheme,  like 
that  suggested  for  the  RQR  [35],  means  a  regular  thinned  ar¬ 
ray  which  still  preserves  the  entire  aperture.  Hence,  some 
grating  lobes  will  appear  in  the  subarray  pattern  although  the 
condition  number  of  the  subarray  phase  matrix  is  strongly  re¬ 
duced.  The  effect  of  this  undersampling  can  be  interpreted  as 
a  penetration  of  nonsuppressable  clutter  power  into  the  sys¬ 
tem  through  these  grating  lobes.  From  the  theory  of  subarray 
design  [38, 39],  it  is  known  that  the  use  of  a  random  subarray 
configuration  can  reduce  these  periodicities  due  to  the  multi¬ 
plicative  superposition  of  the  ambiguities. 

For  a  sufficiently  large  array,  it  is  possible  to  find  a  trans¬ 
formation  by  randomly  re-sorting  the  columns  of  the  identity 
matrix.  For  example  E  €  {0, 1  is  the  randomly  gener¬ 
ated  permutation  matrix  which  can  be  divided  into  nonover¬ 
lapping  block  matrices  E  =  [Ei ,  E2, . . .]  with  E,-  G  {0, 1 
and 


EfE,  =  I 


Im  i  =  j 

0  else. 


The  transformation  matrix  might  then  be  determined  via  a 
sum  of  a  small  number  of  these  blocks  T  =  E,-.  Figs. 

12  and  14  show  the  resulting  SNCR  of  the  MTP  for  the  case 
in  which  the  transformation  was  generated  as  a  sum  of  three 
blocks,i.e.  T  =  Ei  -f  E2  +  E3.  It  can  be  seen  that  the  perfor¬ 
mance  of  MTP  is  slightly  worse  than  EVP  for  the  asymptotic 
as  well  as  the  finite  sample  size  case.  Moreover  almost  no 
difference  can  be  recognized  for  different  dimensions  of  the 
clutter  subspace.  Figs.  13  and  15.  The  number  of  operations 
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to  estimate  the  interference  subspace 

R  =  X  (X^Ei  +  X^E2  +  X^Es) 

is  now  NKM  complex  multiplications,  therefore  requiring  only 
2KM  more  complex  additions  than  the  YBP  or  RQR  methods. 
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Matrix  Transformation  based  Projection,  SS  =  12,  K  =  24 
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Figure  12 :  Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for  MTP, 
dimension  of  subspace  12. 


Matrix  Transformation  based  Projection,  SS  =  12,  K  =  24 


Figure  13 :  Histogram  of  SNCR  for  MTP  at  normalized  frequency  0.25,  calculated 
over  10000  trials,  dimension  of  subspace  12. 
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Histogram  ^  Normalized  SNCR  [dB] 
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Matrix  Transformation  based  Projection,  SS  =  24,  K  =  24 
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14 :  Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for  MTP, 
dimension  of  subspace  24. 


Matrix  Transformation  based  Projection,  SS  =  24,  K  =  24 


Figure  15 :  Histogram  of  SNCR  for  MTP  at  normalized  frequency  0.25,  calculated 
over  10000  trials,  dimension  of  subspace  24. 
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3.3.4  Weighted  Projection  -  Lean  Matrix  Inversion  (LMI) 

This  section  describes  a  simple  yet  computationally  highly  efficient  ap¬ 
proach  for  adaptive  spatial  interference  suppression  which  was  origi¬ 
nally  developed  as  a  robust  method  for  use  in  the  presence  of  dispersive 
channel  errors  [40].  The  method,  called  lean  matrix  inversion  (LMI), 
is  based  on  a  weighted  projection.  The  corresponding  interference  sub¬ 
space  is  estimated  by  the  MTP-technique. 

It  is  known  that  the  dispersive  channel  errors  included  in  data  measured 
with  an  array  of  sensors  will  cause  a  decorrelation  between  the  chan¬ 
nels  which  leads  to  a  loss  of  the  low-rank  property  of  the  corresponding 
covariance  matrix  [41].  The  resulting  spread  of  the  eigenvalue  distri¬ 
bution  can  be  interpreted  as  a  leakage  of  interference  power  from  the 
jammer  subspace  into  the  noise  sub-space.  The  clear  separation  of  the 
two  subspaces  disappears  as  a  result  of  this  leakage.  Since  clutter  eigen¬ 
values  have  a  similar  distribution,  the  application  of  this  approach  to  the 
MTI-problem  appear  to  hold  promise. 


The  problem  of  choosing  the  correct  dimension  of  the  jammer  subspace 
occurs  as  a  result  of  the  eigenvalue  spread  [41]  (see  also  section  3.3.1). 
As  mentioned  in  section  3.3.1,  it  is  obvious  that  the  clutter  power  con¬ 
tained  in  the  remaining  eigenvalues  will  not  be  suppressed  when  the 
dimension  of  the  clutter  subspace  is  chosen  to  be  Af  =  12.  Asa  result, 
a  SNCR-loss  occurs  over  the  entire  frequency  range.  It  has  been  noted, 
however,  that  the  loss 


P 


N 


(28) 


is  strongly  dependent  on  the  clutter-to-noise  ratio  and  probably  on  the 
orientation  of  the  array,  i.e.  side-looking  versus  forward-looking  system. 
The  latter  is  the  subject  of  future  work.  Let  the  inverse  of  the  space-time 
covariance  matrix  in  eq.  (28)  be  substituted  by  its  eigen-factorization 
R“‘  =  UA“'  U^,  where  the  unitary  eigenvector  matrix  is  separated  into 

one  part  S  G  G  corresponding  to  the  eigenvalues  As  G  R  and 

N-x-iN-M)  ,  .  (N-M)x{N-M) 

a  part  N  G  G  corresponding  to  the  rest  An  G  R 

Then,  the  SNCR  in  the  sidelobe  region  of  the  clutter  is  given  as 


(SAg ‘S"  -b  NA^’N")  d  d"  ^NA^’N")  d 

p=  _  = 


N 


(29) 


For  An  close  to  the  identity  matrix,  the  SNCR  would  be  p  =  1  outside 
the  main  clutter  with  S^d  =  0.  With  the  eigenvalue  distribution  shown 
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in  Fig.  2  it  is  evident  that  strong  clutter  with  large  additional  eigenvalues 
(i.e.  far  beyond  the  noise  level)  results  in  an  additional  loss  to  the  SNCR. 
Or,  differently  formulated,  relatively  weak  CNR  produces  eigenvalues 
which  are  at  or  below  the  noise  floor  level. 

The  subspace  dimension  could  be  increased  to  ensure  that  all  dominant 
eigenvalues  are  covered,  however,  this  can  lead  to  a  loss  in  antenna  gain 
as  the  remaining  degrees  of  freedom  in  the  array  may  be  too  small  to  al¬ 
low  detection  of  the  moving  targets.  In  other  words,  a  larger  clutter  sub¬ 
space  dimension  results  in  an  unwanted  suppression  of  target  power  it¬ 
self.  Hence  one  may  conclude  that  no  subspace  dimension  can  be  found 
for  which  projection  methods  work  as  well  as  those  methods  based  on 
the  Wiener  solution. 

The  loaded  sample  matrix  from  eq.  13  can  be  approximated  as 

R  +  yI=  S(As+yI)S^  +  N(An+YI)N" 

^  S(As  +  yI)S^  +  YNN^ 

=  s(As-fyI)s"+y(i-ss") 

=  SAsS^+Yl-  (30) 

Using  the  matrix  inversion  lemma,  and  the  fact  that  the  eigenvectors  are 
mutually  orthonormal,  the  inverse  of  the  loaded  space-time  covariance 
matrix  can  be  written  as 

(R-fyl)"' ^I-SGS^.  (31) 

MxM 

The  entries  of  the  diagonal  matrix  G  =  diagjGi  i , . . . ,  Gmm]  ^  IR 
are  given  as 

Gmm  ~  "b y) 

Eq.  (31)  shows  that  the  LSMI  approximates  a  "weighted”  projection, 
wherein  the  "importance”  of  each  subspace  base  vector  (eigenvector)  is 
determined  by  the  magnitude  of  its  corresponding  eigenvalue.  This  ap¬ 
proximation  motivates  a  modification  of  the  classical  projection  meth¬ 
ods  in  order  to  make  them  more  robust  against  an  eigenvalue  spread. 

As  shown  in  section  3.3.3  the  QR-decomposition  of  the  estimated  sub¬ 
space 

R  =  X  (X^E,  -F  X^E2  +  X^Es)  =  QL  =  Q[QsQn]  (32) 
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can  be  used  to  both,  orthonomialize  the  subspace  base  vectors  and  to  get 
an  approximation  of  the  jammer  eigenvalues  needed  for  the  weighting, 
i.e.  the  diagonal  elements  Fig-  16,  for  instance,  shows 

die  distribution  of  the  magnitudes  of  the  diagonal  elements  of  L  in  eq. 
(32).  These  magnitudes  are  a  measure  of  the  mutual  linear  independence 
of  the  columns  in  R.  The  number  of  dominant  diagonal  values  is  the 
same  as  the  eigenvalues,  i.e.  24.  The  special  shape  reflects  the  inherent 
structure  of  the  asymptotic  covariance  matrix  and  might  be  useful  in  the 
determination  of  the  clutter  subspace. 

The  weighting  matrix  can  be  determined  as  diag(G;w»  =  Lm/{Lm+y))- 
Only  those  diagonal  elements  of  the  QR-decomposition  that  are  larger 
than  a  threshold  a  are  used.  The  threshold  is  the  criterion  to  determine 
the  dimension  of  the  clutter  subspace.  The  value  can  be  chosen  in  iden¬ 
tical  fashion  to  LSMI  as  a  =  3a^,  where  describes  the  receiver  noise 
power.  This  method  was  termed  lean  matrix  inversion  (LMI)  method. 
To  what  extent  the  LMI  is  sensitive  to  the  choice  of  the  threshold  a 
should  be  examined  in  future  work. 


Figure  16 :  Magnitude  of  the  diagonal  elements  after  QR-decomposition  of  the 
space-time  covariance  matrix. 

Figs.  17  to  20  show  the  loss  of  SNCR  for  the  two  different  clutter  sub- 
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space  dimensions,  12  and  24.  It  can  be  seen  that  a  dimension  of  12 
is  not  large  enough  to  show  the  impact  of  the  weighting.  The  results 
in  that  case  are  similar  to  the  MTP;  increasing  the  subspace  dimen¬ 
sion  to  24  allows  the  utilization  of  the  information  contained  in  the  tail 
and  shows  superior  performance  compared  to  all  other  techniques.  It 
should  be  pointed  out  that  the  SNCR-loss  in  that  case  is  even  less  than 
with  LSMI.  This  is  remarkable  since  the  performance  of  LSMI  and  EVP 
were  believed  to  represent  the  optimum  for  the  application  of  interfer¬ 
ence  suppression.  To  verify  that  this  effect  is  not  due  to  variations  in  the 
estimation  of  p,  the  number  of  trials  in  Figs.  18  and  20  was  increased  to 
100000,  thereby  obtaining  a  better  estimate  and  a  smoother  histogram. 


To  investigate  the  performance  of  the  LMI  on  the  form  of  the  transfor¬ 
mation  matrix,  the  three  block  matrices,  T  =  £?_j  E;,  have  been  ran¬ 
domly  generated  in  Fig.  21  for  one  fixed  set  of  estimated  data  vectors 
(snapshots).  For  this  particular  set  of  snapshots,  with  an  SNCR-loss 
of  -1.57dB,  the  variance  due  to  the  random  transformation  matrix  is 
about  lOdB  smaller  than  the  fluctuation  of  the  data  itself.  This  result 
shows  that  the  suppression  performance  of  LMI  and  the  MTP  are  not 
very  sensitive  on  the  choice  of  T. 


Figure  17:  Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for  LMI, 
dimension  of  subspace  12. 
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Normalized  SNCR  [dB]  ^  Histogram 
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18 :  Histogram  of  SNCR  for  LMI  at  normalized  frequency  0,25,  calculated 
over  10000  trials,  dimension  of  subspace  12, 


Lean  Matrix  Inversion,  SS  =  12,  K  *  24 


Lean  Matrix  Inversion,  SS  »  24,  K  »  24 


Normalized  Frequency 


Figure  19 :  Signal-to-clutter  plus  noise  ratio  versus  normalized  frequency  for  LMI, 


dimension  of  subspace  24. 
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-10  -9  -8  -7  -6  -5  -4  -3  -2  -1  0 

Normalized  SNCR  [dB] 

Figure  20 :  Histogram  of  SNCR  for  LMI  at  normalized  frequency  0.25,  calculated 
over  100000  trials,  dimension  of  subspace  24. 

Lean  Matrix  Inversion,  SS  =  24,  K  =  24 


Figure  21 :  Histogram  of  SNCR  for  LMI  at  normalized  frequency  0.25,  calculated 
over  10000  trials,  dimension  of  subspace  24.  Applied  on  one  fixed  esti¬ 
mated  covariance  matrix  but  with  varying  transformation  matrices. 
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4.  Computational  Complexity 


This  section  presents  compares  the  computational  complexity  of  the  different 
techniques  examined  in  this  report.  Only  complex  operations  (multiplications 
and  divisions)  are  counted  for  each  technique  as  they  represent  the  most  intensive 
calculations.  Their  number  is  dependent  on  the 

•  number  of  channels:  N, 

•  number  of  pulses;  M, 

•  number  of  samples;  K 

•  and  the  dimension  of  clutter  subspace:  C. 

To  determine  the  computational  overhead  of  the  SMI/LSMI-method,  one  must 
first  consider  the  estimation  of  the  space-time  covariance  matrix  R  =  ^ 

+  51.  Only  the  upper  (or  lower)  triangle  has  to  be  calculated  as  Ae  co- 
variance  matrix  is  hermitian.  The  number  of  complex  multiplications  is  then 
KNM{NM+\)  ^  Instead  of  a  direct  inversion  of  R  it  is  numerically  more  economi¬ 
cal  to  solve  the  equation  Rvv  =  d,  e.  g.  via  a  Cholesky-decomposition.  Since  R 
is  hermitian  positive  definite,  there  exists  a  unique  lower  triangular  L  with  pos¬ 
itive  diagonal  entries  such  that  R  =  LL^.  Since  the  matrix  L  is  triangular,  the 
solution  can  be  found  by  a  two-step  forward-backward  substitution  procedure 
ljb  =  d  and  L^w  =  b.  The  number  of  complex  operations  required  to  factor¬ 
ize  the  covariance  matrix  is  The  two  substitutions  give  an  additional 

—  NM  operations.  Consequently,  the  overall  number  of  operations  is  of 
the  order  0(^  +  ^). 

A  full  eigenvector  factorization  for  the  EVP  technique  needs  approximatively 
0{\ON^M^)  complex  operations  [32]. 

The  estimation  of  the  clutter  subspace  S  within  the  MTP  or  LMI  technique  is 
achieved  by 


S  =  X  (X"Ei  +  X"E2  +  X^Ea) .  (33) 

This  operation  requires  NMKC  complex  multiplications  since  the  permutation 
matrices,  E,,  only  select  the  appropriate  colunms  of  X^,  they  do  not  calculate 
them.  The  orthonormalisation  of  S  then  requires  another  NMC^  operations. 

Inserting  sample  values  for  the  simulation  (i.e.  N  —  12,M  =  12,X  =  24  and 
either  C  =  12  or  C  =  24)  yields  for  the  number  of  complex  operations: 
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•  SMI/LSMI;  746496, 

•  EVP:  29.87 10^ 

•  MTP/LMI:  4 1 472  if  C  =  1 2  and  82944  if  C  =  24. 

For  this  scenario  the  MTP/LMI  techniques  are  about  10  —  20  times  faster  than 
the  LSMI,  and  about  375  -  750  times  faster  than  the  EVP. 
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5.  Conclusions 


In  this  report,  two  different  recently  developed  fast  rank  reduction  techniques,  the 
so  called  subspace  or  projection  methods,  have  been  applied  to  the  problem  of 
clutter  suppression  for  MTI.  The  basic  theory  behind  them  was  briefly  reviewed 
and  their  small  sample  size  performance  was  studied  via  simulations. 

As  a  result  of  this  study  the  following  important  conclusions  can  be  drawn: 

•  Both  subspace  techniques,  namely  MTP  and  LMI,  are  capable  of  suppressing 
the  clutter  effectively.  The  loss  of  SNCR  between  the  LSMI  and  EVP  on  one 
hand  and  the  MTP  and  LMI  on  the  other  hand  is  only  a  fraction  of  one  dB. 

•  It  is  not  necessary  to  suppress  the  full  clutter  subspace  defined  by  Brennan’s 
rule,  particularly  in  the  finite  sample  case.  It  is  the  magnitude  of  the  eigenval¬ 
ues  rather  than  the  number  larger  than  the  noise  floor  which  is  the  dominant 
factor,  i.e.  their  contribution  to  the  overall  clutter  power.  This  reduction  leads 
to  savings  in  computational  overhead. 

•  The  asymptotic  and  finite  sample  size  behaviour  of  each  method  are  notice¬ 
ably  different.  Ripples  in  the  shapes  of  the  matched  filter  outputs  for  known 
covariance  matrices,  blur  or  even  disappear  in  the  different  realizations  due 
to  the  uncertainty  of  the  estimation  process. 

•  The  LMI  technique  was  superior  to  all  other  methods,  in  that  it  possessed 
both  the  smallest  loss  of  SNCR  and  the  lowest  computational  complexity. 
The  first  point  is  especially  noteworthy,  since  LSMI  and  EVP  were  believed 
to  represent  the  optimum  interference  suppression  performance. 

Since  this  study  is  based  only  on  one  particular  simulation  scenario,  future  work 
will  have  to  focus  on  different  open  items,  such  as: 

•  Investigation  and  comparison  of  the  proposed  techniques  for  varying  number 
of  channels  and  pulses;  particularly  for  the  practical  case  of  a  small  number 
of  channels  with  a  larger  number  of  pulses. 

•  Studying  the  impact  of  a  different  radar  orientation,  such  as  a  forward-looking 
system  (i.e.  range-dependency  of  the  clutter),  as  well  as  the  effects  of  inher¬ 
ent  channel  errors  on  the  performance  of  MTP  and  LMI,  especially  for  small 
sample  size. 

•  Developing  a  method  for  automatic  determination  of  an  optimum  dimension 
of  the  clutter  subspace  for  EVP,  MTP  and  LMI,  e.g.  based  on  the  ratios  of 
the  eigenvalue  magnitudes. 
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•  Comparison  of  the  performance  versus  computational  cost  of  the  fully  adap¬ 
tive  processor  concept  with  a  so-called  beamspace  approach.  In  the  latter  the 
data  are  dimension-reduced  before  calculating  the  cancellation  filter.  This 
technique  has  also  been  proposed  to  reduce  the  computational  effort  [1]. 
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Annex  A 

Proof  of  Eg.  (20)  _ 

NxM 

Let  Q  €  (D  be  an  unitary  matrix  and  A  a  diagonal  matrix  with  real  entries 
>  ■  ■  ■ )  then: 


maxtrjQ^AQ}  = 


(A.1) 


Proof: 

The  trace  in  eq.  (A.l)  can  be  re-written  as 

N  N 

tr{Q^AQ}  =  tr{AQQ"}  =  £ 

m=l  m=l 

where  ^  denotes  the  m-th  column  of  the  identity  matrix  In-  Since  QQ^  de¬ 
scribes  a  projection  matrix  (i.e.  it  is  idempotent),  it  has  M-eigenvalues  equal  to 
1  and  the  rest  0.  Hence  it  yields  for  the  diagonal  elements  am  = 

i)  0  <  ttm  <  1  m=l,...,N 


ii)  £a;;,=  tr{QQ"}=M. 


Proposition: 


N  M  N  M 

£  Xmam  =  £  XmOmi  +  £  XmOm  <  £  Xm 

m=l  m=l  m=:Af+l  m=l 

or 

M  N 

£  Xm{am  —  1)  ■+•  £  XmOm  ^  0, 

m=l  m=A/+l 

respectively.  Since  A-m+i  is  greater  than  or  equal  to  zero,  it  follows 

M  N 

£  (^m  ~  ~  1)  +  £  <  0  (A.2) 

m=l  m=^M+\ 


M  N  /  M  N 

XmiS^m  1)  XmO-m  ^  X\f-\-\  I  (O-m  1)  "b 

m=l  m=A/+l  \m=l 

=  ^M+l  f  £  Otm  -  • 


(A.3) 


Using  ii)  the  term  within  the  brackets  on  the  right  side  of  inequality  (A.3)  equals  zero. 

(q.e.d) 
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